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LOCATING THE DISCONTINUITIES OF A BOUNDED FUNCTION 
BY THE PARTIAL SUMS OF ITS FOURIER SERIES I: PERIODICAL 

CASE 


GEORGE KVERNADZE, THOMAS HAGSTROM, AND HENRY SHAPIRO 


Abstract. A key step for some methods deeding with the reconstruction of a function 
with jump discontinuities is the accurate approximation of the jumps and their loca- 
tions. Various methods have been suggested in the literature to obtain this valuable 
information. 

In the present paper, we develop an algorithm based on identities which determine the 
jumps of a 2r-periodic bounded not-too-highly oscillating function by the partial sums of 
its differentiated Fourier series. The algorithm enables one to approximate the locations 
of discontinuities and the magnitudes of jumps of a bounded function. We study the 
accuracy of approximation and establish asymptotic expansions for the approximations 
of a 2x-periodic piecewise smooth function with one discontinuity. By an appropriate 
linear combination, obtained via derivatives of different order, we significantly improve 
the accuracy. Next, we use Richardson’s extrapolation method to enhance the accuracy 
even more. For a function with multiple discontinuities we establish simple formulae 
which “eliminate” all discontinuities of the function but one. Then we treat the function 
as if it had one singularity following the method described above. 


1. Introduction 

It is well known that the main difficulty in applying a Fourier series as a tool for approx- 
imating a discontinuous function is the Gibbs phenomenon. Namely, the approximation 
of a function by the n-th partial sum of its Fourier series is only of order 0(1 jn) for each 
point of continuity of the function and oscillations are 0(1) in an 0(l/n) neighborhood 
of the discontinuity point. 

Two distinct approaches to resolve this difficulty have been suggested in the literature. 
The first is to reduce the oscillatory behavior by filtering. The second is to use step func- 
tions to reconstruct the discontinuous function. The latter approach was first suggested 
by Gottlieb et. al. [17] and has been further developed in [1], [2], [5], and [16]. The key 
step in the method of reconstruction suggested in [5] is the accurate approximation of 
the location and the jumps of a given function. 

1991 Mathematics Subject Classification, 65D99, 65T99, 41A60. 

Key words and phrases. Locating discontinuities, Fourier series, asymptotic expansions. 

The second author was supported in part by NSF Grants DMS-9304406 and DMS-9600146. 


1 



2 


GEORGE KVERNADZE, THOMAS HAGSTROM, AND HENRY SHAPIRO 


Later, Eckhoff [10], [11] considered a different approach to locate the discontinuities 
using Prony’s method. As a result he developed an efficient method of approximating 
the locations of singularities and the jumps of a piecewise smooth function with multiple 
discontinuities. The approximations are found as the solution of a system of algebraic 
equations. 

To justify the importance of allocating the discontinuities and the jumps of a function, 
let us give a brief review of the idea of reconstruction of a function from its truncated 
Fourier series as developed in the above mentioned papers. 

Let g be a 27r-periodic function which is piecewise smooth on the period with a finite 
number, M, of jump discontinuities. In addition, we assume that the first 2n + 1 Fourier 
coefficients of the function are known. If G{6 ) = (ir — 0)/2, 8 £ (0,27r), is the 27r-periodic 
sawtooth function, then the assumption that the function g is piecewise smooth on [— 7P, it] 
with a finite number of singularities is equivalent to the following representation of the 
function: 

i M—l 

(1) «(«) = - E (9]-»G(« - «»)+»(*). 

* m=0 

where 8 m and [<7]m> m = 0,1,. .. ,M — 1, axe the locations of discontinuities and the 
associated jumps of the function g , and g is a 27r-periodic continuous function, which is 
piecewise smooth on [— 7r,7r]. 

Hence, the problem is to find a good approximation for the constants 8 m and [</] m , 
given the first 2n + 1 Fourier coefficients of the function g. Then g can be recovered from 
the partial sums of its Fourier series based on identity (1) and the undesirable Gibbs 
phenomenon could be avoided. 

Recently another approach to recovery a piecewise smooth function was suggested by 
Geer and Banerjee. (See [4], [13], and [14].) The authors introduced a family of periodic 
functions with “built-in” discontinuities to reconstruct a piecewise smooth function with 
exponential accuracy. The main assumption of the method is knowledge of the jumps and 
the locations of discontinuity of the given function. To find these, the authors suggested 
the following: use the well-known formula of symmetric difference of the partial sums 
of Fourier series which determines the jumps of a bounded function to obtain the first 
estimate for the location of discontinuities; then utilize the modified least-squares method 
to improve the accuracy of approximation. It should be mentioned that a method for 
the recovery of a piecewise smooth function with exponential accuracy, utilizing the 
Gegenbauer polynomials, was developed in a series of papers by Gottlieb and Shu (see 
[18] and the indicated references). But again, the authors assume some knowledge of the 
location of the singularities of the function. 

In the present paper, we consider an essentially different approach for the approxima- 
tion of the points of discontinuity and the jumps of a function based on special formulae 
determining the jumps of a bounded not-too-highly oscillating function by the partial 
sums of its differentiated Fourier series. It is shown that the largest local maximum of 
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the absolute value of the differentiated partial sums of the Fourier series occur in the 
vicinity of the actual points of discontinuity of the function. Furthermore, for a piece- 
wise smooth function with one jump discontinuity, we establish asymptotic expansions 
for the approximations of the location of the discontinuity and the magnitude of the 
j um p Utilizing the expansion formulae, we use Richardson’s extrapolation method to 
achieve higher accuracy. For a function with multiple singularities, we establish simple 
formulae which “eliminate” all discontinuities of the function but one. Then we treat the 
modified function as if it had only one discontinuity, using the method described above. 


2. Definitions 


Throughout this paper we use the following general notations: N, Z+, Z, and R are the 
sets of positive integers, nonnegative integers, integers, and real numbers, respectively. 
L[a, b] is the space of integrable functions. W[a,b] is the space of functions on [a,b] 
which may have discontinuities only of the first kind and are normalized by the condition 
g(8) = ( 0 ( 0 +) + g(8—))/2, 8 E (a, b). (Here, and elsewhere, ff(0+) and g{8 — ) mean the 
right and left hand-side limits of a function g at a point 8, respectively). C[a,b] is the 
space of continuous functions on [a, 6] with uniform norm || • ||[ a ,6]- By C p [a, b], p E N, 
we denote the space of p-times continuously differentiable functions on [a, 6]. 

All functions below are assumed to be 27r-periodic with the obvious exceptions. 

If g E L[— 7 r, 7 r], then g has a Fourier series with respect to the trigonometric system 
{ljcosnfljsinn#}^^, and we denote the n-th partial sum of the Fourier series of g by 
Sn{g,0), i-e., 

S n (g ; 8) = + ]T( a fc(s) cos k8 + b k (g ) sin k8), 

1 k = 1 


where 


If* If*. 

a k (g) = — / g(r) cos krdr and b k (g) = — / ^(t) sin krdr 

7T J — ir 7T J-n 

are the fc-th Fourier coefficients of the function g. 

By S n (g\ 8) we denote the n-th partial sum of the conjugate series, i.e., 


Sn{g\8) = Y^i a k{g) sinks - b k (g) cos k8). 

k=i 


Correspondingly, by g we denote the conjugate function, i.e., 

m = Hm{— 1 r g(a + ;)-» l l ? . z l W 

cr 7 t Jh 2 tan t 


which exists and is finite almost everywhere for any g E L[— 7r, vr] (cf. [19, Theorem, p. 
79]). 

By K we denote positive constants, possibly depending on some fixed parameters and 
in general distinct in different formulae. For positive quantities A n and B n , possibly 
depending on some other variables as well, we write A n = o( B n ), A n = 0(B n ), or 
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A n ~ B n , if limn—oo A n JB n = 0, sup n€jV A n /B n < oo, or < A n /B n < K 2 , respectively, 
where K\ > 0 and K 2 > 0 are some absolute constants. 

Definition. Let A = (A k)T=i be a nondecreasing sequence of positive numbers such 
that 


( 2 ) 


OO 1 

k = i 


A function / is said to have A-bounded variation on [a, 6 ], i.e., / € ABV[a,b], if 


sup£ 


n 


k = 1 


!/(»>) - /(**)! 

Ait 


< oo, 


where II is an arbitrary system of disjoint intervals ( 0 ^, 6 *,) € [a, b). 

We say that a function / is of harmonic bounded variation on [a, 6 ], i.e., / € iffJV'fa, 6 ], 
if A* = k, k € JV\ 

Remark 1 . For a reader unfamiliar with Af?V[a, 6 ] classes of functions we give some 
basic properties of these classes. 

The A- variation “measures” the total oscillation of a bounded function. ABV[a,b] 
is a generalization of Vfa, b], the class of functions of bounded variation (obviously 
A BV[a, b] = V[a, &] if A fe = 1, k € N). 

Waterman [23, p. 108] mentioned that the inclusion 


(3) ABV[a,b] C W[a,b] 

holds for any A BV[a,b] class of functions. 

It is known as well [22, Theorem 3, p. 114] that for rBV[a,i>] and A BV[a,b] Wa- 
terman’s classes of functions, defined by the sequences T = ('yk)kL\ a^d A = (A^^Lj, 
respectively, the inclusion ABV[a, 6 ] C r5V[a, i] holds if and only if ££=1 I/ 7 k = 

1/A*). 

The constraint on the sequence A is natural, since if series ( 2 ) converges, ABV [a, 6 ] = 
B[a , b ] , where B[a, b) is the class of all bounded functions on [a, b]. This makes it clear that 
the HBV[a,b] class is sufficiently wide and “almost” covers B[a,b], since A = (A: 1 +£ )£L 1 
converges for any e > 0 . □ 

If there is no ambiguity, we shall usually suppress the dependence on the domain and 
simply write C, ABV, etc. 


3. Main identities 

The identity determining the jumps of a function of bounded variation by means of 
the partial sums of its Fourier series has been known for a long time: 
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Theorem 1 ([8] and [12]). Let g £ V. Then the identity 
(4) lim^ki^) = I (9 ( #+ )_ s (<l_)) 

V J n— oo n 7T 

is valid for each fixed 6 E [ — tt, 7r] . 

Golubov [15] generalized identity (4) for Wiener’s [24] V p classes of functions and higher 
derivatives of the partial sums of Fourier and conjugate series. Further generalizations, 
extending the results of Golubov to A BV classes of functions, have been obtained by one 
of the authors. 


Theorem 2 ([20]). Letr E Z+ and suppose KBV is the class of functions of K-bounded 
variation determined by the sequence A = Then 

a) the identity 


(5) 


n-+oo n^ r+1 


(-ix 

(2 r + l)7r 


(,(*+) -,(«-)) 


is valid for every g £ KBV and each fixed 6 £ [ — 7r , 7r] if and only if 


( 6 ) 


KBV C HBV. 


b) There is no way to determine the jump at the point 6 £ [ — tt, 7r] of an arbitrary function 
g £ AJBV by means of the sequence [S^ T \g] £)) • 


Theorem 3 ([20]). Let r £ N and suppose KBV is the class of functions of K-bounded 
variation determined by the sequence A = Then 

a) the identity 


(7) 


lim 

n— ► oo 


5l !r) (s; *) 


71 


2 r 


(-l) r+1 

2r7r 


(*(«+)-*(«-)) 


is valid for every g £ KBV and each fixed 6 £ [— 7T,7r] if and only if condition (6) holds, 
b) There is no way to determine the jump at the point 8 £ [— 7r,7r] of an arbitrary function 

g £ KBV by means of the sequence 0)) . 


Remark 2. Theorems 2 and 3 (see [20, Theorems 1 and 4]) implicitly include the 
following statement: if g £ C fl HBV, then the convergence of (5) and (7) to zero is 
uniform with respect to 9 £ [ — 7r, 7 t] . □ 

Furthermore, as a simple corollary from Theorems 2 and 3 follow the identities which 
determine the jumps of the derivatives of a continuous function. 
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Corollary 1. Let r E N and r — p be a positive odd number, and suppose g E C p 1 is 
such that g ^ E HBV . Then the identity 


( 8 ) 


lim 

n— i-oo 


sSHsi j) 

n'-r 




is valid for each fixed 6 E [— 7r,7r]. 

Proof. By virtue of (3), g ^ E HBV C W. Hence 

( 9 ) = 

for r > p. Then identity (8) instantly follows from (9) and Theorem 2. □ 
The following statement is proved similarly. 


Corollary 2. Let r € N and r — p be a positive even number, and suppose g E C p 1 is 
such that g&) E HBV . Then the identity 


(10) lim = i 2 V p >(0+) - g (p \9-)) 

v y n r ~ p (r-p) 7T 

is valid for each fixed 9 E [— 7r,7r]. 


4. Preliminaries 


In what follows we need the following additional notations. 

By 9 m = B m {g) and [g) m = g{9 m +) - g{9 m ~), m = 0,1,. . . ,M - 1, we denote the 
points of discontinuity and the associated jumps of a function g E W. By M = M(g) we 
denote the number of discontinuities (finite or infinite) of the function g E W. 

For a fbced p E Z + , r E N, and g E L we set 


( 11 ) 


DT n (p\r]g\6) 


(r — p) 7T 

n r ~P 


(— 1 ) L -%~ L Sfr\g; 6) ifr — p is odd, 
(— l) I i £-1 5^(5; 6) if r — p is even. 


For a fixed p E Z + , r E N, and M E N, the points 9 m {p\T\ g\n), m = 0, 1, . . . , M — 1, 
axe defined via the following condition: 


(12) \DT n (p; r; g\ 0 m (p; r; g\ n))| = max{|DT n (p;r;p;0)| : 9 E B{9 m \ A{g))} , 

where B(9 m ) A(g)) is the closed bail around 9 m with the radius A(g) = |min{|0 m — 
9k | mod 2i r : m, k = 0, 1 , . . . , M — 1 and m ^ fc}. 

To simplify notations, we sometimes omit fixed parameters and write DT n {9), DT n {g\ 9), 
or DT n (r\g\9). Similarly we simplify the notation for tf m (p;r;p;n). 

By G(9) = (tt — 9)1 2, 9 E (0,27t), we denote the 27r-periodic extension of the sawtooth 
function. If 7 E R, then following the notations in [11] we set 

G(r, 9) = G{9 - 7) and G* +1 ( 7 ; 9) = J G k (r, 0)d6 
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for Jfc E Z + , where Go = G and the constants of integration are successively determined 
by the condition 

[ Gk{T,T)d,T = 0 . 

It is trivial to check that 


(13) 


S'n(G-, 9) 



Dn(9) 


1 

2’ 


where 

(14) 


D n {9) = \ + 5 Z COS ke = 
1 k=\ 


for* S' 2 

n + | for 5 € 27rZ 


is the Dirichlet kernel. 


Lemma 1. Let r £ N be fixed. Then 

a) the closed form of the following sum exists: 


(15) 


71 




rh"' tI+ k + 5"” ,+ 


5 


where the last term contains either n or n 2 . 

b) The following expansion holds for every a — 1 € N : 


(16) 


n— 1 

£ 


k = 1 


1 

F+i 


= C(r + 1) 


l i_ _ l l 
rn T 2n r+1 


+ E(-i) 


5-1 


B s T(r + 5) 1 


+ 0 (- 


s= 2 


(2^)! T(r + l)n r+s n 


,r+a+l 


). 


where £(r) = Hfclj k~ T , r > 1, is the Riemman zeta function, T is the Gamma function 
and B s , s € N, are Bernoulli numbers. 


Statement a) of the lemma can be found in [9, p. 1]. Using the Laplace method, 
the proof of expansion (16) is a simple corollary of the integral representation of the 
Hurwitz zeta function [3, Theorem 12.2, p. 251] and Watson’s lemma [6, p. 253]. It was 
generously offered by Prof. E. Coutsias [7]. 


Lemma 2 (Bernstein’s inequality). IfT n is a trigonometric polynomial of degree n 
N, then 


ll T nll[a,6] < 7 ||T n || [a) i], 


6 


where [a, b] C [— 7r,7r]. 
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Lemma 3. Let a function g € C q be such that € V. Then 

a) g 6 C 9-1 and g ^~ 1 '> e Lip a for all a € (0, 1), i.e., |fir^ 9_1 ^(^) — ^^ _1 )(r)| < K \8 — rl 0 
for some K > 0 and all 6 ,t £ R. 

b) The following estimates hold: 

(17) Ms), Ms) = o(^), 

where Rn{g) = ||5 , „(^; •) - g\\[- T ,sr) and Rn{g) = ||<S n (0; •) - g\\[-n,r], n £ N. 

Proof. Statement a) can be found in [19, exercise 3, p. 81]. As regards statement b), 
by virtue of Holder’s inequality, since g € C q , we have: 

Ms), Ms) < E(Mg)l + Ng)l) = f: Mg( ' ))l , H ; 

k=n k=n K 

( oo -i \ !/ 2 / oo \ V 2 

E^J (E(a fc (^) 2 + M^) 2 )j . 

Meanwhile, it is known [21] that if g € C H V, then 

00 1 

(is) EWj ) 2 + M?) 2 ) = °(-)' 

k = n 71 

Now (17) follows as a simple combination of (16), (18), and (19). □ 

The following axe some basic properties of the function D^(8), r G N. 

Lemma 4. Let ip n = <p n ( r ) > 0 and ij) n = ipn( r ) > 0 be the closest nonzero roots to the 
point zero of the equations D^ T \d) = 0 and D^ r+1 ^(8) = 0, respectively. Then for any 
fixed r G Z+ : 

*)***&%)• 

c){-iy+'D(™\<p n )~n™. 

d) (— l) r+1 I?^ 2r+1 l(^) is increasing on [— ip n {r + l),y> n (r + 1)], concave on [— < p n (r + 1), 0] 
and convex on [0,<^„(r + 1)]. 

e) ( — 1 ) T D^ r \8) is a 2 ~k - periodic even and smooth function with the global maximum 
attained at 8 = 0. In addition , the sequence of the absolute values of the local maxima is 
decreasing as a function of 8 € [0,7r] and 

(20) |D< 2 '>(0)l > JT(r)|D!r>(,i„)|, 

where K(r) > 1 and n 6 N. 
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Proof, a) Let us prove the statement for an even n, i.e., n = 2 n. By (14) we have 


S ignD<r(£) = sign ((-1)' |£> cos ± k*’cos£ 


fn - 1 

E 

it=i 


2n 


fc=rn+l 


2n 


/ /"-I fc x n ~ 1 

= sign ^(-l) r ^ fc2r cos ^ + E( 2n ~ cos 


(2n — fc)7r 


2n 


= sign ( (-l) r (j2(k 2T - (2 n - A:) 2r ) cos ^ - (2n) 2r 

Vfc=i 


k-K 


2 n 


( 21 ) 


= sign( — l) r+1 . 


Again by (14), sign^n^^) = sign(— l) r for n € N and 0 6 [0,7r/4n]. The latter com- 
bined with (21) and the Mean Value Theorem instantly guarantees <p„ € (rjAn, 7r/2n). 
Similarly we treat the case when n is odd. 

b) The statement is proved analogously and we omit the details. 

c) According to (14) and (15) 


(22) (-i y + 1 Dg r+ 1 ) ( 8 ) = 2 fc2T+1 sin *0 < E fc2r+1 - n2r+2 - 

k = i fc=i 

Meanwhile, since <p n € [7r/2n, ir/n] (see a)), taking into account the well-known inequality 
20/7T < sin 0 < 0 for 0 G [0,7r/2], we have 


n o [ n /2] 1 [ n /^ l 

(23) £ ifc 2r+1 sin(Jfc<p«) > - £ * 2r+ Vn > 1 E fc2r+2 - " 2r+2 > 

kml * km 1 n fc=l 

where [a] means the integer part of a number a. Combination of (22) and (23) completes 
the proof of statemant c). 

d) Since the function (-1 ) r+1 D£ 2r+2 )(0) is positive on [-<p n (r + l),<p n (r + 1)] (see (14)), 
(— l) r+1 £)£ 2r+1 )(0) is monotonic on the interval. Furthermore, (— l)( r+1 ^£)^ 2r+3 ^(0) is pos- 
itive and negative on [-rp n {r + 1), 0] and [0,tp n (r + 1)], respectively. But (p n (r + 1) < 
V'n^+l) (see a) and b)). Hence (-l) r+1 £>) l 2r+1) (0) is concave and convex on [-ip n {r+l), 0] 
and [0 ,<p n (r + 1)], respectively. 

e) Let us prove inequality (20) as the rest of the statement is trivial. It is clear that 


(24) 



for n > 4. But by virtue of (15), lim n _*oo q n > 1 as well. The last combined with (24) 
implies the existence of K(r) such that 


(25) q n > K(r) > 1 

for n > 4. 
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Figure 1. n = 16. 

= |[ E + E ) k 2T cos ktp n \ 

\^V / n€[^/2,3?r/2] kiffn£[ir/2,3ir/2]J 
n 

< E * 2r > 

k>n/4 

since ip n satisfies the estimate b) and the sums in (26) have different signs. The rest 
instantly follows from (24)-(26), and the identity 0) = (— l) r Hk=i ^ 2r - Validity of 
(20) for n < 4 is trivial. □ 

5. General idea of algorithms and the accuracy of approximations 

This is the general idea of all the following algorithms: according to identities (5) and 
(7), if g E HBV, then for a fixed r £ N , p = 0, and sufficiently large to € JV, the 
function \DT n (6)\, 0 G [— 7r,7r], (see (11)) must attain the largest local maximum nearby 
the actual points of discontinuity of the function g, since at the the points of countinuity 
of g, DT n {6) = o(l) by virtue of Theorems 2 and 3. (The proof of Theorem 4 includes 
a rigorous proof of this statement.) Hence we search for the singularity locations of a 
function by locating the largest local spikes of the differentiated partial sums of its Fourier 
series. 

Figures 1-4 represent the graphs of the normalized differentiated partial sums ^S n (g] •) 
of the function (65) with increasing n. They illustrate the dynamics of creation of sharp 
spikes in the vicinity of the actual points of discontinuities of the function. 


Besides 

|ff>Wn)l 

( 26 ) 
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Figure 3. n = 64. 
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Figure 4. n = 128. 

Now we study how well the points 8 m (n) and the values DT n (8 m (n )) approximate the 
points of discontinuity 8 m and the jumps [g] m of a function g. 

5.1. Approximation to the points of discontinuity. Let us first consider the worst 
possible case. 


Theorem 4. Let p = 0 and r 6 N be fixed, and suppose g £ HBV is a function with a 
finite number, M, of discontinuities. Then the estimate 


(27) 


8 m (n) = 8 m + 


1 


o(i) 


[$]mA(s) n 


is valid for each fixed m = 0, 1, . . . , M — 1. 


Proof. Without loss of generality let us make several assumptions. We assume that 
M = 2 and r = 2r+ 1 is an odd number. The points of jump discontinuity of the function 
g are 8 0 = 0 and 8\. We shall prove estimate (27) for 8 0 as it is completely analagous for 
8\ by virtue of the periodicity of g. 

Now let us set 

(28) g{8) = g{8) - [ -^-G(8) - ^G{8 X - 8). 

7 T 7T 

It is obvious that 


(29) g e C n HBV, 

since continuity of g follows from (28). Moreover, since G € V C HBV and HBV is a 
linear vector space (see [23, p. 108]), g € HBV as well. 
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Besides by virtue of (11), (13), and (28) 

DT n { g] 0 ) = ^DT n (G ] 0)+ [ ^DT n (G(0 l] -) ] 9) + DT n (g ] 0) 

_ ( — l) r (2r + l)[g]o (2f) 

n 2r +i n ^ ' 

+ ( ~ ir(2 2 r ^ 1 — D<? r \d - 6 0 + DT n {g\ 6) 


(30) = Io(n] 6) Ii(n\ 9) ER(n\0). 

It is obvious that |/ 0 (n; 6 ) | attains the global maximum at 6 = 0 and without 7i(n; 6) 
and ER(n\ 9) terms we could exactly locate the discontinuity point Oq = 0 just searching 
for the global maximum of \DT n (g\ 0)| on the period. By virtue of (29) and Remark 2, 
ER(n\ 0) contributes a small error independent of 0 £ [— tt, 7r], i.e., ER(n\ 0) = o(l). But 
according to (14) and (30) 

(-l) r (2r + l)[g]q f sin((n + |)(g-gi)) y 2r) 

2 sin J 


h{n\9) = 


n 


2r+l 


V 


( 1),( „t« 1)bll (E <£(» + \t ■*.((» + 1)(« - «.) + y ) 

. / , l > 2r sin (( 7l +|X 0 - 5 l)+ r7r ) 

+ ( n + o ) 


k - 0 
(2r-Jfe) 


[5] 


,2 sin 2 

^(r) 


2 sin 


0-0, 


) 


(31) ^ A ( ,rv 

as well for 0 £ B(0; A(g)). Hence 

(32) e n = ||7i(n; •) + ER{n\ OII[-A(»),A(y)] = °(!)- 
Consequently, by virtue of statement e) of Lemma 4 and (32), we have 


(33) |Jo(n; 0)| - e„ > |7 0 (n; <p n ) | + e n 

for sufficiently large n £ N. But (33) combined with (12), (30), and statements a) and 
e) of Lemma 4 already guarantees 

|*o(n)| < <p n <- 
n 

for sufficiently large n £ N. 

Next, to achieve a more accurate estimate, namely (27), we use a simple estimate of a 
root of an equation. 

First, let us mention that since 0o(rc) is the extremum point, then 

(34) DT' n (g ; 0 o (n)) = 0, 
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which itself implies (see (30)) 

(35) /o(n; 0 o (n)) = 0 o (n)) - ER'(n\ 0 o (n)) = T n (6 0 {n)), 

where T n is am n-th degree trigonometric polynomial. 

According to estimate (32), (35), and Lemma 2 we have 

(36) ll T n||[-A(s),A( 9 )] = 

Let us assume for simplicity that [^jo > 0. Furthermore, we know /g(n; 6) is odd 
decreasing and convex on [— ip n (r + 1 ), 0] and concave on [0,y> n (r + 1)]. (See statement d) 
of Lemma 4.) Hence the line passing through the points (±y> n (r + 1), / 0 (n; ±<£> n (r + 1))) 
will occur below the positive part of the function J o (n;0) and above its negative part. 
So, for sufficiently large n £ N, 6o(n) will satisfy the inequality 

(37) |0o(n)| < |So(n)|, 


where 0o(rc) is the solution of the following equation 


(38) 


7p(n; y n (r + 1)) 

<Pn(r + 1 ) 


e = T n {6). 


Here the left hand side of the equation represents the above mentioned line. 

Hence, by virtue of (11), (13), (30), (36), and statements a) and c) of Lemma 4, we 
obtain 

* o(n) - 

which combined with (37) completes the proof. □ 

Let us now consider a more typical case. 


Theorem 5. Let p = 0 and r £ N be fixed, and suppose g is a piecewise continuous 
Junction such that g £ HBV. In addition, we assume that M(g) and M(g ) are finite. 
Then the estimate 

(39) M*n) = «.(,) + ^ (l9'l~0(l) + S® El*^)) 

is valid for each m = 0,1,... , M(g ) — 1. 

Proof Again for simplicity let us assume that M(g) = 2, M(g ) = 1, r = 2r + 1 is 
odd, and # 0 ( 5 ) = # 0(5 ) = 0. Furthermore, let us introduce a function g now via the 
following identity: 

m = 9(6) - - - i|s'loG,(«). 

X m=0 7r 


( 40 ) 
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Since the conditions of Theorem 5 in particular imply the conditions of Theorem 4, by 
similar arguments we conclude that 6o(n) satisfies the estimate (37), where $o( n ) is the 
solution of equation (38) now with 

T.(t) = -/,(»;>)-( 

(41) = -^DTUO(e l] )-,e)-^DT‘ n (G l -,e)-DT‘ n (g;«)- 
Hence, to complete the proof it is enough to estimate T„. 

By the construction g € C fl HBV (see (29)). Consequently, according to (11), (41), 
and Remark 2 

(42) II-E-® (ra; *) II [-*■>*] ~ °(^)* 

The estimate for I\ directly follows from (31). Namely, 

(43) lllt-iCIAM = ^jO(l). 

As regards I§\ by virtue of (11), (13), (15), and (41) we have 

(44) ll(4 1) )'(n; Ollt-^l = [5]o0(l). 

The combination of (37), (38), (41)-(43), and (44) completes the proof. □ 

Now we turn our efforts to study probably the most interesting case: a 27r-periodic 
piecewise smooth function with one jump discontinuity. As expected, the approximation 
in this case is significantly more regular. Namely, the following statement holds. 


Theorem 6. Let p — 0 and r £ N be fixed, and suppose the function g piecewise belongs 
to C\ q> 2, and has a single discontinuity at 8o G (— 7r,7r). In addition , we assume that 
g(q) £ y. Then there exist constants K\, A2 ,. . . >K q such that 


(45) 
Namely 

(46) 


K x K- 


K a 


6q{t\ n) = 6 q + — + — + . . . ^7 + °(^r)- 


n 


n 


n 


n 


9+1 


r + 2 [g'] 0 
71 id] 0 


and K 2 


r 2 [g ]q 

r lff]o 


for r > 2. 

In particular , if the derivative of the function g does not have a jump at 0 q, then the 
approximation has the order 0(1 /n 4 ). 
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Proof. Let us first assume that r > q. We shall establish an algorithm for computing 
the constants K \ , K 2 , • • • , K q , and perform the actual computation for K\ and K 2 . 

Without loss of generality we assume that 60 = 0, r = 2r + 1, and q = 2q + 1. Now we 
consider the function g defined by 

1 2?+l 

(47) m = «(«) - - £ (9 (i) ]oGr t («). 

* teo 

Since the function g in particular satisfies the conditions of Theorem 5, by virtue of 
(39) there exists a constant Kq such that 

(48) l«o(n)| < ^ 
for n £ N . 

As we know (see (34) and (47)), 6o{n) satisfies the following identity 

1 + l 

(49) DT' n (g-, 0„(n)) = - £ «o(n)) + DT^g; 0„(n)) = 0. 

* k= 0 

By construction ^( 29+1 ) £ C D V C C fl HBV. Hence by Remark 2 and Lemma 2 we 
have 

(50) si 2r+2) (g; 0) = S£2r+l-2 ? )(-(2, + l). ^ = ^Jr+1-2,) 
uniformly with respect to 6 € [ — tt, 7r] . 

Furthermore, expanding expression (49) into a Taylor series around zero on the interval 
[— K 0 /n 2 , K 0 /n 2 ] and taking into consideration (11), (13), (48), and (50), we obtain: 

bloPf «>(O)0„(n) + ii)< J '«)(0)«oW 3 + |4 Sr+6, (0)9c,(n) s + . . . 

+ l D^ 2 q+ 2 r+ 3 ) (a 0 ) ( 2Jr °) 2<?+2 ) 

+ (2 g + 2)! n n 4 “+ 4 } 

+is']o(i>r > (o)+iDr j >(o)«„w + + • ■ • 

+ (2 9 + l)! n { ^' n) n 4q+2 J 

+ b (!,+1) ]o(4>‘ 2r - !,, (0) + 

(51) + o(n 2r+1_2? ) = 0, 
where \nk,n\ < K 0 /n 2 , fc = 0,1, . . . ,2q + 1. 

It follows from (15) that all error terms in the Taylor expansion have order 0(n 2^-2, ). 
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The expression for I>i r) (0), r E Z+, (see (15)) suggests seeking an expression of 6 0 (n) 
in the form (45). 

According to equation (51), since the error term has an order o(ra 2r ~ 2?+1 ), all coefficients 
of n k , k > 2r — 2q + 1, must equal to 0. This condition generates the set of equations 
with respect to the yet unknown constants K\,K 2 , ■ • ■ , K' 2 q+i- 

We set up one by one the equations for powers of n, with decreasing order of de- 
gree, starting from n 2T+1 . It is clear that by (15) and (51), only two terms, namely 
[0]oZ?i 2T+2 ^(O)0o(7i) and contribute n 2r+1 and n 2r . Consequently, the com- 

parision of the coefficient leads to the following system of linear equations with respect 
to Ki and K 2 (see (14), (15), and (45)): 


(-l) r+1 [s]o 


n 2r+3 Ki 
2r + 3 n 2 


+ (-l)V]o 


n 


2r+l 


2r+ 1 


= 0 


and 


(-l) r+1 [<7]o 


Ki' 


' n 2r+3 Ki n 2r+2 

2r + 3 n 3 2 n 2 


n 


2 r 


+ (~i) r [$ ] 0_ 2~ — 


which instantly implies (46). 

Furthermore, let us observe that the highest degree of n contributed by each term of the 
sequence Q t = (l>l 2r+2 '' 2,) (0)^o(n) 2 ' _ ‘ _1 ) 2 io 1 , Z = 1, 2, . . . , q + 1, ignoring the constants 
of expansion, is 2r — 21 + 3. 

Now we proceed by induction. Let us assume that the constants K \ , if 2 , • • ■ ,K 2 i-3, 
and K 2 i- 2 are already defined by setting up equations with respect to the coefficients of 
n degree less then 2r — 2(Z — 1) -f 3. Next, we shall show that a new system of equations 
for the coefficients of n 2r ~ 2l+3 and n 2r ~ 2l+2 represents a system of linear equations with 
respect to K 2 i- 1 and K 2 \. In addition, the determinant of the system is nonzero , and 
hence the system is consistent. 

Indeed, the only terms which may contribute K 2 i-i and K 2 i unknowns axe in the 
sequences Qj, j < l. Hence, by (15) and (45) we have 


£) (2r +2 i— 2 i)( O )0 o ( n )2i— * 


n 


2r+2j-2t+l 


+ lower degree terms 
2r + 2j — 2i + 1 j 


(52) 


(K, Kv-x , 

x ( — + ... + — — + 


Vn 2 


n 


21 


k 2 , 

n 2l+l 


+ 


°(^>) 


2 jr — t — 1 


Consequently, the highest degree of n contributed by this product with factor K 2 i-\ is 

,2r+2j— 2t+l /if 1 \ 2 J-‘- 2 K 2 l-i 




2r + 2 j -2 i + 


i(§y 


n 


2 ( 


1 ~ n 2r-2/+3+2(l-j) 


But, 2r — 21 + 3 + 2(1 — j) < 2r — 21 + 3 unless j = 1. Hence only the sequence 
Qj = (Z)| l 2r+2 ^(0)^o(n), D( l 2r '(0)} contributes the constant K 2 i- 1 and it clearly appears 
in the first degree in the expression for 8o(n) . (We treat the case for K 2 i similarly.) In 
addition, the determinant of the linear system with respect to K 2 i - 1 and K 2 i is triangular 
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with nonzero diagonal entries, (— l) r+1 [<7]o/(2r + 3) 7^ 0 and (— l) r+1 [5)0/2 ^ 0, and that 
guarantees the solvability of the system. 

Finally, the equation for n 2r-2i+3 defines K21-1 , so the equation for n 2r-2,+1 will 
define iC^+i- Let us mention that the coefficients Ki,K 2,... , .K29+1 depend only on 
Wo,b']o [j (2s+1, ]o. 

To prove the theorem for the case r < q we need some minor changes in the arguments. 
First, one has to use expansion (16) for D| l 2r+2 -' _2, ^(0) in estimate (52) as soon as 2r + 
2 j — 2i < 0. Second, to estimate the error term in (49), i.e., (50), let us mention the 
following: since y^ 2r+2 ^ G C 2q ~ 2r ~ 1 , then by virtue of (17), expanding g into Taylor series 
around zero on the interval [— Ko/n 2 , Ko/n 2 ], we have: 


sr + 2 ) (s; *o(n)) = + 

(53) 

which, ignoring the constant factor, represents the desired estimate for the error term. 
The rest of the proof is completely analogous. □ 

Taking an opportunity of the lucky similarity between the coefficients (46), using a 
simple linear combination of expansion (45) for f and r, 2 < f < r, we significantly 
improve the accuracy of approximation. Namely, the following statement holds. 


Corollary 3. Let p = 0 and suppose a function g piecewise belong to C q , q > 3, and has 
a single discontinuity at 60 G (— 7r,7r). In addition, we assume that g(i) £ y Then for 
each fixed f and r £ N , 2 <f < r, there exist constants K\, K2 , . . . ,K q such that 


(54) 


r(r + 2) . . r(r + 2) . K\ 

2(r^ o(rin) ' 2F^) o( ;n) = 9o + ^ + 


+ 






5.2. Approximation to the jumps. Now let us study the approximation to the jumps 
of a function. 


Theorem 7. Let p — 0 and r G N be fixed, and suppose g is a piecewise continuous 
function such that g G HBV. In addition, we assume that M(g ) and M(g ) are finite. 
Then the estimate 

DT n ( g] e m (n)) = [g} m + 0(-) 

n 

is valid for each m = 0, 1, . . . , M(g ) — 1. 
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Proof. Again for simplicity we assume that M(g) = 2, M(g') = 1, r = 2r + 1 is an 
odd number, and # 0 ( 5 ) = 9o{g) = 0. By virtue of (40) we have 

DT n ( g] d 0 (n )) = I^PT n (G;5o(n)) + ^r n (G(^i-)^o(n)) 


7 T 


(55) 


bio 


7 r 


OT n (G i; flo(n)) + Dr B (j;flo(»)) 


Further analysis is trivial as we take the Taylor expansion of (55) around zero on the 
interval [— K 0 /n 2 , Ko/n 2 ]. By virtue of (11) and (13) we get 

OT n ( S ; «„(»)) = (Df)(0) + 

+ ( ~ 1),( n t« 1)M ' °£ irl w n )-^ 1 ) 

+ ( ~ 1),(2 I,i l ' 1 - ^ 1 ° Pf- 11 (».(„)) 


(56) 


« (-l) r (2r + 1 )tt (2rV , 

+ ^+i 5 n (5^o(n)), 


where |i/„| < K 0 /n 2 . 

Taking into account that g 6 C fl HBV, the rest of the proof follows from (14), (15), 
(31), (56), and Remark 2. □ 

Now, am interested reader will easily fill out the details of proof for the following 
statement. 

Theorem 8. Let p = 0 and r G N be fixed, and suppose a function g peicewise belong to 
C\ q > 2, and has a single discontinuity at 9 q G (— 7T,7r). In addition, we assume that 
g(<i) £ y . Then there exist constants K\, K 2 , ■ ■ ■ , K q such that 


(57) 
Namely 

(58) 
and 

for r > 2. 


DU r; *o(r;n)) = b]o + ^ + § + ■•• + % + <K^)- 


n rr 


= j bio 


„ r 2 r + 2 , r + 2 bio 

= ijblo + — b lo - — 


Extrapolating expansion (57) in r based on identity (58), we improve the accuracy of 
aproximation. Namely, the following statement holds. 
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Corollary 4. Let p = 0 and suppose a function g piecewise belong to C q , q > 2, and has 
a single discontinuity at 6 0 G (— 7t,7t). In addition, we assume that g^ G V. Then for 
each fixed f and r G N, 2 < f < r, there exist constants K\, Ki, . . . ,K q such that 

(89) «o(n)) - -f—.DUr, «»(-*)) = bio +$ + ••■ + ^ + <K^)- 

r — r r — r vr n 9 n 9 

5.3. Approximation to the discontinuities and the jumps of derivatives of a 
function. As a simple corollary of (8), (9), (10), and Theorems 6 and 8 we obtain 
estimates for the location of points of discontinuity of the derivatives of a continuous 
function and the associated jumps. Below we represent just typical statements. 

Theorem 9. Let a function g G C p ~ l , p G N, piecewise belong to C q , p + 2 < q, and 
suppose g( p> has a single discontinuity 6 0 G (— 7r,7r). In addition, we assume that g* q ' > G V. 
Then for each fixed r G N there exist constants K\ , K 2 , ■ ■ . , K q ~ v such that 


(60) 

0o (p;r;g;n) = 6 0 (g (p) ) + 

*■ + 
7l 2 

1 Kq-P . / 1 \ 

” n q - p+1 K n q - p+ 

Namely, 





_ r - p + 2 [s (p+1) ] 0 
1 r-p [gWjo 

and 

r - p + 2 [$ (p+1) ] 0 
! r-p b«]„ 

for r — p > 2. 





Theorem 10. Let a function g G C v 1 , p G N , piecewise belong to C q , p + 2 < q, 
and suppose g^ has a single discontinuity 8 0 G (— 7 r, 7 r). In addition, we assume that 


e v . 

Then for each fixed r G N, there exist costants K\,Ki, 

» . . j Cq~.p such ihai 


DT n (p-, r; g ; 6 0 (p\ r; g\ n )) = [s (p) ] 0 + + . . . + ? V 

n n 9 p 

+ o( ). 

Namely , 




Jfi = r -^\s ir) ]o 


for r — p 

> 2. 



6. Description of the algorithm 

Now we describe the main idea of the algorithm which we propose to locate the points 
of discontinuity. For simplicity, we assume that the function is piecewise smooth. 

In case the function has a single discontinuity, according to Corollary 3 we search for the 
global maximum of \DT n {6)\ for fixed f and r,2 <r < r. Afterwards, utilizing expansion 
(54) and applying Richardson’s method of extrapolation we improve the accuracy. 

The situation drastically changes if the function has more than one point of disconti- 
nuity. In this case we do not have expansion (54) for the approximation. 
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To overcome this difficulty we generate, for a fixed r (E N, the sequence of partial sums 
of Fourier series of functions (fif m )m=o > defined via the recursion relation 

(61) g m + 1 ( 0 ) = (1 - cos(0 - 9 0 {g m \n))) d g m {8), 

where g 0 = g, d E N is fixed, and 9 G {g m \n ) is the location of \DT n (g m ] 0)|’s global 
maximum on the period for sufficiently large n € N. 

The idea behind generating this sequence is to “eliminate” (at least numerically) the 
discontinuities of the function g one by one by simply searching for the global maximum 
of the function \DT n (g m ; 0)| on the period. In other words, “removing” the highest jump 
of the function, we can see the second largest jump. 

A straightforward computation based on simple trigonometric identities generates the 
Fourier coefficients of g m +i via the Fourier coefficients of g m . Below we represent the 
identities for d = 3 (see (61)): 

= 2 

15 15 

- — cos^o(Sm;«)(o/t-i(firm) + Ofc+i(^m)) + y sin 8 0 (g m ; n){b k -i (g m ) - b k+1 (g m )) 

3 3 

+ - cos26 0 (g m ]n)(a k . 2 (gm) + ak+2{9m)) - ^ sin20o(Sm;n)(&fc-2(5m) - b k+2 {9m)) 

- i COS 39o(g m -,n)(ak-3(gm) + O-k+sigm)) + £ sin30 O (5m;^)(^-3(5m) - ^+3^)), 

o o 


— 

15 15 

- — cos0o($m;n)(5*_i(0 m ) + 5 fc+ i(ffm)) - — s'm9 0 (g m -,n)(a k -i(g m ) - a k+1 (g m )) 

O O 

3 3 

+ -cos 20 o(g m \n)(b k - 2 ( 9 m) + h+ 2 (g m )) + -sm28 0 (g m -,n)(a k . 2 (g m ) ~ a kJr2 (g m )) 

- ^ cos3do(Pm;n)(6fc_ 3 (5m) + h+3 (9m)) ~ ^ sin 28 0 (g m ] n)(a k - 3 (gm) - a k+3 (g m )), 

O O 

where a k (g m+1 ) = a- k (g m+ i) and b k (g m+1 ) = -&_ fc ( 5 m+ i), k £ N. 

Our approach to eliminate the points of discontinuity could be justified by the following 
observation: multiplying a function by the factor (1 — cos (# — 9o(g m i n ))) 3 ' w e are not 
adding a new point of discontinuity but significantly reducing the jump of the function 
g m at 6 m . More precisely, if 6 m — 8(g m ]n) = 0(n~ p ), then the function g m +i (see (61)) 
at the point 8 m will have the jump of (1 — cos(0 m - 8 0 (g m \ n))) d [Sm+i]o ^ n~ 2dp [g m+ i ] 0 . 

Figures 5-7 illustrate a step by step “elimination” of the points of discontinuity of 
the function (65) utilizing formulae (62) and (63). They correspond to the graphs of 
^S'n(9m]-) for m = 1,2,3. 


®fc(9m+ 1 ) 


( 62 ) 

and 

bk{9m+ 1 ) 


( 63 ) 
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Figure 5. n = 128. 



Figure 6. n = 128. 
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Figure 7. n = 128. 

Summarizing all the above we suggest the following algorithm, which we have imple- 
mented using Mathematica: 

Steps 1-3 find initial approximations for the discontinuities of the function g. These are 
then refined in Steps 4-5. 

Step 1 (Initialization): 

Select some fixed d € N, which will be used throughout the entire algorithm. 

Let r = 1; this value is used in defining DT n (g m ] 0) in Steps 2 and 3. Let 2 + 1 
coefficients of the function g be given. Select some small value for n 0 (usually 
no = 16); this value should be a power of 2 and is used as the starting subscript 
of the sequence constructed in Steps 2a and 2b. Let m = 0 and go = g. 

Step 2 (Find a discontinuity, if one exists): 

Step 2a (Find the maximum of |Z>T„(g m ;0)| over the period): 

Step 2a(l): Using the adaptive plotting routine internal to Mathe- 
matica, with the number of initial points set at the maximum of 25 
(the default) and (to ensure that we do not miss the maximum), 
determine the point constructed by the plotting algorithm which has 
the largest value of \DT no (g rn \ 0)|. 

Step 2a(2): Using the 0-value of the point found in Step 2a(l), ap- 
ply Newton’s method to find the 0-value where the maximum of 
\DT no (g m \6)\ occurs as the solution of the equation DT no (g m ; 0) = 0. 

Let us denote the 0- and correspondinglDT^^m; 0)|-value by 0(n o ) 
and r(no). 

Step 2b (Determine if the 0- value found in Step 2a is a discontinuity): 
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Form the sequences 0(n o ), 0(2no), 0(4no), . . . and r(no), t(2tio), r(4no), 
... by successively doubling n 0 and finding the maximum of \DT n (g m ] 0)|. 
Find the maximum by using Newton’s method, using the previous 0-value 
in the sequence as the initial approximation. 

Stop when one of three conditions arises. 

Condition 1: The ratio of two successive r-values is less than 0.6 
(see (57)). The point is not a point of discontinuity. Go to Step 4. 
Condition 2: The estimated relative error between two successive 
0- values is less than some predefined tolerance. The point is a point 
of discontinuity. Go to Step 3. 

Condition 3: The number of terms in DT n (g m \ 0) exceeds 2ni. The 
point is a point of discontinuity. Go to Step 3. 

Step 3 (Remove the discontinuity): 

Let 0 m denote the final value of 0 determined in Step 2b. Increase m by 1, 
define the Fourier coefficients of the function <7 m +i(0) = <? m (0)(l — cos (0 
utilizing (62) and (63), and return to Step 2. 

Step 4 (Refine the estimates of the points of discontinuity): 

Two conditions can arise: no discontinuities were found in Steps 2 and 3 (termi- 
nate the algorithm) or one or more discontinuities were found (continue). Let M 
be the number of located discontinuities and 0o, 0i, ... , 0 a/-i their locations. 

Select some 2 <t < r and some n > no, n = 2 P , to be used in Step 4a. 


Do the following for m = 0 to M — 1. 

Using the values of 0 m current at the time form the function 
9m( 6 ) = 9(0) n;=o]#m (! - cos (0 - 0i)Y> and apply the extrapolation 
method of Step 4a. 

Step 4a: 

Form the sequences 0(f;no), 0(f;2no), 0(f;4no), • • ■ , 0(r;n) and 0(r;no), 

0(r; 2no), 0(r ; 4no), . . . , 0(r; n), which represent the 0-values where \DT n (g^\ 0)| 
takes on its maximum value. This can be done efficiently by using New- 
ton’s method with 0 m as the starting point and solving the equation 

DTfa„ ;t)=0. 


Using the extrapolation defined by equation (54), form the sequence 

K f + 2) r(r + 2 )^_, 

0( n ) o7I — ~ oTI — ^( r > n ) 


2(r — f) 

for n = n 0 , 2no, 4n 0 , ... , n. 


2(r — r) 


Perform Richardson’s extrapolation on the sequence 0(no), 0(2no), 0(4no), 
. . . , 0(n). 

Replace 0 m with the final value obtained. 
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Step 5 (Additional refinement of the estimates): 

If only one discontinuity was detected in Steps 1-3, then stop. Otherwise repeat 
Step 4 a second time. 

Our implementation uses Mathematica’s capability to do high-precision arithmetic 
and track the loss of precision that occurs due to accumulating roundoff error. Besides 
displaying the final answers it displays the number of digits of precision that remain. 

7. Some numerical results 

In order to illustrate the numerical results obtained by the described algorithm, we 
will consider several examples. 

First we consider a 27r-periodic function g with two discontinuities such that g € HBV , 
namely 

n(S) =1° if-7T<0<O 

^ ' ( 6 sin |+2 ifO<0<7r 

It is obvious that the function g is not piecewise absolutely continuous and therefore 
the methods suggested in [4], [5], and [11] fail. Still according to Theorem 2 we should 
obtain an o(l/n) approximation. We find the absolute value of the largest error for 
approximation to the points of discontinuity as follows: 


N 

32 

64 

128 

256 

Location- error 

EBB3I 

2.2(— 2) 

6.4(-3) 

2.0(— 3) 


The following piecewise smooth function was considered in [5]. 


(64) 


9(8) = 


sin | if 0 < 6 < 0.9 
— sin | if 0.9 < 9 < 2tt 


Below we present a detailed description of all computations. The first table shows the 
error in the approximation to the location of the discontinuity by differentiated Fourier 
partial sums of degree r = 7 and r = 8, and then their linear combination via formula 
(54). 



r = 7 

r = 8 

linear combination by (54) 

2 

2.42( — 1) 

2.40( — 1) 

1.82(— 1) 

4 


iaRiTBHl 

1.90(— 2) 

8 


1.81(— 2) 

1.50(— 3) 

16 

4.90(— 3) 

4.76(— 3) 

1.06(— 4) 



1.22(— 3) 

7.05(— 6) 

64 

3. 19( — 4) 

BEDES 

4.54(— 7) 

128 

IQ 

IQ* 

O 

00 

co" 

00 

bo 

00 

T 

3 
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Here we present a full table of Richardson’s extrapolation started from the last column 
of previous table. 


N 

r = 7,8 







2 

l-S(-l) 








1.9(— 2) 

8.0(— 3) 






8 


3.4(-4) 

8.9(— 5) 





16 

l-0( 4) 

1.2(-5) 

2.3(— 6) 

9.2(— 7) 




32 

7.0(— 6) 

4.4(-7) 

4.6(— 8) 

1.0(— 8) 




64 

4-5(— 7 

1.4(-8) 

8.2(— 10) 

9.6( — 11) 

1.2(-11) 

1.8(— 12) 


128 

2.8(— 8) 

4.7(— 10) 

1.3(— 11) 

7.7(— 13) 


EEBDB 

2.1( — 14) 


The following is the error in the approximation to the jumps of the function using 
r = 7 and r = 8, and their combination (59). 



r = 7 

r = 8 

linear combination by (59) 

2 

2.86(0) 

3.35(0) 

5.60(— 1) 



1.22(0) 

1.74{-1) 

8 


5.17(-1) 

4.55(— 2) 

16 

2.20(— 1) 

2.37(— 1) 

1.14(— 2) 

32 

9.90(— 2) 

1.13(— 1) 

2.87(— 3) 

64 

4.85(— 2) 

5.55(— 2) 

-4 

l — 1 

OO 

I 

128 

2.40(— 2) 

2.74(— 2) 

1.79(— 4) 


This table is Richardson’s extrapolation applied to data above. 


N 

r = 7, 8 







2 

5.6(-l) 







4 


4.5(-2) 






8 


2.6(— 3) 

3.5(— 3) 





16 

l.l(-2) 

1.2(-4) 

to 

to 

T 

l-K-5) 




32 

to 

bo 

1 

CO 

2.9(— 6) 

l-4(— 5) 

4.1(-7) 

61(— 8) 



64 

7.1(— 4) 

BSffil 

9.3(— 7) 

l-l(-8) 

17(-9) 



128 

l-7(— 4) 

l-OC-7) 

5.8(— 8) 

3.0(— 10) 

5.6(— 11) 

2.9(— 11) 

8.3(— 12) 
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The next example has been considered in [11]. 


(65) 



0 

e* 

cos | 
0 


if 0 < 6 < 1, 
if 1 < 6 < 2, 
if 2 < 6 < 5, 
if 5 < 6 < 2ir 


Applying the suggested algorithm the following results have been obtained. Here we 
simply present the largest error in approximating any of the three discontinuities. 


N 

32 

64 

128 

256 

512 

Location-error 

1.7(— 4) 

6.1C-7) 

1-4C-8) 

3.5( — 11) 

9.7(— 14) 


8. Conclusion 

Let us give some comments on our results. 

As we already mentioned, the formula which determines the jumps of a bounded not- 
too-highly oscillating function by means of its Fourier series has been known for a long 
time. But to our best knowledge it has never been utilized for a numerical approximation 
of the locations of discontinuity points. 

Theorems 2 and 3 state that it is possible to detect the locations of discontinuities 
and the jumps under the condition that a bounded function does not have too high 
total oscillation (condition (6)). On the other hand, Fourier series fail to distinguish a 
continuous functions from discontinuous one, if the function is too highly oscillating (the 
necessity of condition (6)). 

It follows from Theorem 4 that identities (5) and (7) represent a powerful tool for the 
allocation of the points of discontinuity of an almost arbitrary function - excepting the 
minor restriction on the variation of the function and a finite number of discontinuities we 
impose no conditions. Still the approximation is of order o(l/n). The factor ([^] m A(y)) _1 
confirms a logical observation: the smaller the jump of a function and the distance 
between the points of singularity, the more difficult it is to detect its location. 

Taking into consideration asymptotic expansions (45), (54), (57), and (59), we think 
that the method is best suited for a piecewise smooth function with a single discontinuity, 
although applying the suggested method of “removing” discontinuities leads to good 
results for a function with multiple discontinuities too. 

For piecewise smooth functions, we can obtain very high orders of approximation. The 
numerical results confirm that high accuracy is indeed attainable with fairly low degree 
trigonometric polynomials. 

Regarding numerical results, applying expansion (54) for different pairs, we observed 
higher accuracy for larger values of r. For instance, the accuracy of the location of 
discontinuity is only order of 1CT 9 for the function (64) applying (54) and Richardson’s 
extrapolation for r = 2, 3, and n — 128. Numerical results confirmed priority of expansion 
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formula (54) over (45): we gained three digits of accuracy. It should be mentioned as 
well that increasing the order of the partial sums we obtained better results: for the same 
function (64) using expansion (54), r = 6,7, combined with Richardson’s extrapolation 
for n = 512 we achieved an accuracy of 1CT 20 for the location of the point of discontinuity. 

The numerical results were obtained from a program written in Mathematical which is 
available online through the third authors home page http : //www . cs .unm . edu/ shapiro/ 
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